Equivalence of Matrix Models for Complex QCD Dirac 
. Spectra* 
O . 

. G. Akemann 



m 
(N 

(N 
> 



Service de Physique Theorique CEA/Saclay 
Unite associee CNRS/SPM/URA 2306 
F-91191 Gif-sur-Yvette Cedex, France 



Two different matrix models for QCD with a non- vanishing quark chem- 
ical potential are shown to be equivalent by mapping the corresponding 
partition functions. The equivalence holds in the phase with broken chi- 
■ ral symmetry. It is exact in the limit of weak non-Hermiticity, where the 

' chemical potential squared is rescaled with the volume. At strong non- 

, Hermiticity it holds only for small chemical potential. The first model 

' proposed by Stephanov is directly related to QCD and allows to analyze 

■ ■ ' the QCD phase diagram. In the second model suggested by the author all 

microscopic spectral correlation functions of complex Dirac operators can 
'' be calculated in the broken phase. We briefly compare those predictions to 

(~| , complex Dirac eigenvalues from quenched QCD lattice simulations. 
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' 1. Introduction 

Ten years ago the proposal 1 was made to use random matrix models 
in the study of the low energy phase of Quantum Chromodynamics (QCD), 
and it has turned out to be very successful. Since then many new matrix 
model correlation functions have been calculated for that purpose. The field 
theoretic origin for the applicability of matrix models has been understood 
and a wealth of data from QCD lattice simulations has been compared 
to their predictions. We refer the reader to 2' for a review summarizing 
these developments. However, the introduction of a chemical potential /i 
for quarks has left many open questions until very recently, concerning both 
numerical simulations as well as analytical matrix model predictions. The 
reason is that for non- vanishing fi the eigenvalues and thus the determinant 
of the Dirac operator become complex. 



* Presented at the Worksliop on Random Geometry in Krakow May 2003, Poland 
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Let us summarize what can be calculated using random matrix models 
for QCD. The predictions can be put in two groups: the Dirac operator 
spectrum and the phase diagram. In the first case the distribution of Dirac 
eigenvalues Aj rescaled with the volume V, VXi = const., are predicted from 
matrix models close to the origin A = 0. This region is magnified microscop- 
ically as it is most sensitive to spontaneous chiral symmetry breaking. This 
fact can be seen from the Banks-Casher relation, p(0) = V{qq)/'K, which 
relates the Dirac operator spectral density p{X) at the origin to the order 
parameter for chiral symmetry breaking {qq), the chiral condensate. The 
predictions of matrix models are universal and depend only the Nf quark 
flavors, their masses mj, the topological charge v and the global symme- 
try breaking pattern. They can be related to a particular limit of chiral 
perturbation theory, in which only the zero momentum modes of the Gold- 
stone bosons contribute. The matrix model description breaks down at the 
Thouless energy, when non-zero momentum modes come into play. 

In the second application matrix models are used as an effective model 
for QCD to predict the qualitative features of the phase diagram. The order 
and location of the phase transition lines is found by inspecting the non- 
analyticities of the matrix model partition function in terms of parameters 
as the temperature T and chemical potential /x. Here, the close relation 
to field theory is lost and the matrix model resembles a Landau-Ginsburg 
description. The following result was obtained in ^ for QCD with two 
massless quark flavors. A first order (full) line starting at large /x merges 



T 



{qq) = 



(qq) / 




Fig. 1. The schematic QCD phase diagram with Nf — 2 massless flavors from jS]. 

with a second order (dotted) line from p = into a tri-critical point. Other 
transitions from nuclear matter formation or more complicated phases such 
as color super-conductivity have been omitted here, where the latter can be 
addressed using matrix models as well . 
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We take this simplified phase diagram to explain which regions are ex- 
plored. One part of the phase diagram accessible to QCD lattice simulations 
for three colors is the temperature axis at /i = (including more flavors and 
the influence of their mass ratio). The second part is for small values of /i 
in the vicinity of the second order phase transition line up to or crossing the 
tri-critical point. Only in the quenched approximation the full diagram can 
be covered so far. We refer to |S] for reviews and references. On the other 
hand very detailed matrix model predictions for the microscopic Dirac spec- 
trum are available at T = ^ = (see and references therein). For T ^ 
and /i = analytic predictions exist too [HI Ej and have been compared to 
lattice simulations along the T-axis up to and above the transition 8 . 

Here, we will investigate the //-axis at T = instead, for which two 
different matrix models have been suggested [HlEl- The model introduced 
first by Stephanov i9j explicitly includes a chemical potential, using the 
same symmetry arguments as [Q. However, in this model only global prop- 
erties are available, with the microscopic spectrum being unknown. For 
that reason another model ^Uj was proposed where all microscopic spec- 
tral correlators can be calculated. This model is only valid in the broken 
phase by construction. It follows from |lj by analytic continuation of the 
corresponding orthogonal polynomials into the complex plane and contains 
H only implicitly through a non-Hermiticity parameter r. Our aim here is 
to show that ^Uj is equivalent to j9j and thus as closely related to QCD. 
We prove that the two partition functions of these models are equal in 
the broken phase in the limit of infinite matrix size. In the limit of weak 
non-Hermiticity, keeping Vfi'^ fixed, the equivalence is shown for Nj = 1,2 
and 3 mass-degenerate flavors (the case of general Nf and 7^ will be 
treated elsewhere ^3). At strong non-Hermiticity the partition functions 
are mapped for Nf = 1 to leading order in small ji, and we expect the 
equivalence to hold also for A'^^ > 2. The equality of partition functions 
does not automatically imply that all correlation functions are identical as 
well. However, it makes it very plausible that the results for the microscopic 
spectrum |1L)| are common to both models and are universal. We also give 
a brief comparison between the predictions and quenched QCD lattice 
simulations jl21 113j for small /i 7^ in the broken phase. 

2. Equivalence of matrix models for QCD with 7^ 

The following matrix model was introduced by Stephanov ^9^ to study 
QCD with chemical potential 
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Here $ is a complex matrix of size x (N + 1/) with the partition function 
being in the sector of topological charge i^. The quark masses of the 
Nf flavors have to be rescaled with the volume in the same way as the 
eigenvalues in the microscopic limit. The block structure of the Dirac matrix 
follows from chiral symmetry and the location of fi is dictated by adding 
fj,q^q to the quark Lagrangian. The Gaussian weight is chosen such that the 
Banks-Casher relation for /i = is satisfied. In the following we will mostly 
set (qq) = 1 for simplicity. 

The partition function for Nj = 1 quark flavor was already know jl4j . 



Zi{fi- m) = exp[-iV(l + n^)] exp[-iVm^] /q 2Nm Jl + fi' 



(2) 



where Iq denotes the modified Bessel function. Here we have taken already 
the large- limit, assuming that the mass has to be rescaled with an ap- 
propriate power of the volume V = 2N (see eqs. (|15|) . (|2c{|) l. Eq. Q holds 
in the broken phase. At a critical chemical potential fic ~ 0.527 determined 
by the equation 

1 + fil + logfil = 0, (3) 

a first order phase transition takes place. Above fj, > Hc it holds Zf{fj.; m) ~ 
^2N £qj, = q ^2]. This can be show by mapping eq. to a sigma 
model partition function over complex matrices a of size Nf x Nf in eq. 
(|17|) . allowing us to treat also A'^^ > 1 in the broken phase below. The 
macroscopic spectral density of complex eigenvalues z = x + iy given by |H] 



is bounded by the following curve 



(-^2 _ y2)2 



1 



(4) 



4/.2)) + /(x2-4(l-/x2)). (5) 



It splits into two domains only at ^ = 1 > /ic- On the other hand for small 
fi <^ 1 the density simplifies to a constant on an ellipse 



piz-fj.) 



1 



if 



-+4<i- 

4 4//4 - 



Let us compare this to the second model [in]. The quenched Nf 
partition function is defined as 



Zu{t) = J f[(^(fzj w{Zj,Z* 

z*) 



\z\^'+^ew 



N 



(6) 

= 

(7) 

(8) 
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given in terms of complex eigenvalues instead of matrices. The Vander- 
monde determinant is denoted by l^]^{z'^) = Y{k>i{^k ~ ^i)- "^^^ parameter 
r G [0, 1] labels the degree of non-Hermiticity interpolating between real 
and maximally complex eigenvalues for r = 1 and r = 0, respectively. The 
macroscopic density reads^ 



p{z;t) 



27r(l 



if 



+ 



y 



2(1 + t)2 2(1 -t)2 



< 1 . 



(9) 



Identifying the two macroscopic densities eqs. © and @, which implies 
that to leading order in all moments of the two models become equal, 
leads us to^ 



2^^ 



(1-r^) 



(10) 



Clearly the model eq. Q is always in the broken phase as the density 
is constant on a single domain for all r E [0,1). Here, all microscopic k- 
point correlations functions can be explicitly calculated using the technique 
of orthogonal polynomials in the complex plane ^0] (see e.g. eqs. (P7j) . 
H28|) below). For that purpose we define the following complex Laguerre 
polynomials in monic normalization. 



2t\ 



N J 



L 



[ID 



which are orthogonal with respect to the weight eq. 

In order to compare to the partition function eq. above we also have 
to introduce mass terms here, which we do as follows: 



I N Nf 

= (nn(-K 

\i=i /=i 



Zj W(Zj,Zj 



Nf 

/=i 



m}) 




(12) 



Such an expectation value can be calculated using the theorem proved in 
|15j for characteristic polynomials. In our case it only contains orthogonal 
polynomials in monic normalization, and we obtain 



det 



kJ=l,...M 



Zn{T;{mf}) 



N+k-1 



-Nmf 
2t 



/\Nf{m?) 



(13) 



^ Note the missing factor 2 in the corresponding eqs. in |10|. 
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Here, we only display the r-dependent prefactors. The partition function 
defined in eq. (fT^ is thus real as it should. From eq. (fT^ we can take the 
large- iV limit where one has to distinguish now two cases, weak and strong 
non-Hermiticity. 



2.1. The weak non-Hermiticity limit 

In the weak non-Hermiticity limit introduced by ^M. (see also the recent 
review we take the Hermitian limit r 
— > oo at the same time, keeping 



1 and the infinite volume limit 



N{1 - 



a 



(14) 



finite. In this limit the macroscopic density eq. Q becomes localized on 
the real line while the microscopic correlations functions still extend into the 
complex plane. In the microscopic limit we have to rescale the eigenvalues 
and consequently also the masses with the volume. In order to obtain a 
parameter free prediction this is done by 



2N{qq)z = ^ , 



(15) 



and similarly for the masses. Since in the limit r — > 1 the model eq. 
has {qq) = 1/V2 we define the rescaled masses as ^/2Nmf = u>f. Taking 
the asymptotic large- A'^ limit of eq. (|13j) in this way we obtain 



2ii,weakia; Wf}) ~ exp 



det 



(16) 



where we have used that — > exp[— from the weak limit eq. (|14|) . 
The expression eq. (|l()j) only differs from the finite volume partition function 
at /i = by the a-dependent exponential prefactor. 

In order to prove that the same results can be obtained from the matrix 
model partition function eq. we rewrite it in terms of the following 
(T-model representation |14j : 



Zi{^;m) = J da exp[— N (qq)"^ Tr a a'^] det ^ 



a + m 



a^ + m 



N 



(17) 



Here, we have restricted ourselves to mass-degenerate flavors and to the 
sector of topological charge 1^ = 0. Compared to eq. N and Nf have 
interchanged roles, as a is now a complex Nf x Nf matrix. We show ex- 
plicitly for Nf = 2 how to calculate the integral eq. H17|) by parameterizing 
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a. The calculation for Nf = 3 then follows in a similar way. Any complex 
matrix a can be written in the Schur decomposition 



a = U{Z + R)U-^ 



(18) 



where U is unitary, Z =diag(zi, . . . , zat^) is a diagonal matrix containing the 
complex eigenvalues and ii is a complex, strictly upper triangular matrix. 
For Nf = 2 we thus have two complex eigenvalues and R contains a single 
complex number. Setting (qq) = 1 we can write 



(fzi(fz2(fR exp -iV(|zip + |z2|^ + 



I Zip + m(zi + z^) + rr? 



\Z2 - Zl\ 

+ m)R 



N 



X det 



R*{zi + m) 



|z2p + 'm{z2 + Z2) + m? 







driridr2r2dpp I dipidip2d(f) [{r^ - p ){r2 - p ) - p p 



2,,2-\N 



X exp 



-iV(Sj=i,2(r| - 2rjmcos ipj) + 2m^ + p^) 



X [rf + r2 — 2ri r2 (cos ipi cos (/?2 + sin ipi sin ip2 ] 



X exp 



2,,2iAf 



driridr2r2dpp [{r^ - p ){r2 - p ) - p p ] 
-N{rl +rl+ p^+ m^)] [{rl + r|)/o(2A^mri) 



xIo{2Nmr2) - 2rir2li{2Nmri)Ii{2Nmr2)] . 



(19) 



In the first step we have applied the Schur decomposition, which leads to the 
Vandermonde determinant. In the second step we have changed variables, 
Zj Zj = {xj + m) + iyj = rjexp[i(/?j] for j = 1,2, introducing polar 
coordinates here and in i? = pexp[i(p]. In the last step we have carried 
out the angular integrations. So far all steps are exact. Next we employ 
the scaling eq. (fTl)) . p'^ = a'^/{2N). Consequently we can replace the first 
factor in the last equation above by an exponential 

[{rl-p^){rl-p'')-p^p''f ^ {nr2r eM-o^\p^ + rl + rl)/{2rlrl)] , (20) 

using lim7v->oo(l — x/N)^ = exp[— x]. The integral over p thus reads: 

POO 

lim / dpp eyi.v[-Np^][{rl- p^){rl 



2Nrlr. 



2 



exp 



-G?{r\ + rl) 
2r\r2 



(21) 
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If we neglect the subleading o? in the denominator on the right hand side 
the remaining integrals over r\ and in eq. H19|) completely factorize. We 
can thus investigate them separately using a saddle point approximation. 
Because of {qq) = 1 we rescale the mass here as 2Nm = uo. Neglecting 
higher orders in the saddle point is at r?^]^ 2 ~ ^- The two contributions 
containing the factors Iq and Ji can be written in the following form, 



2i,weak{a;uj) ~ exp 



-a 



Io{uj) Ii{uj) 
Ii{uj) Io{u) 



(22) 



This is equivalent to eq. (|16|) for Nf = 2, when taking the limit of degenerate 
masses there. In particular it justifies on the level of partition functions the 
identification eq. (|lflj) , leading to the scaling eq. (|14)) . The same calculation 
can be easily done for Nf = 3, following the same lines as above. In fact, 
recently we managed to show in full generality, that eq. (|17|) leads to eq. 
for all Nf with non-degenerate masses and topological charge v ^ 



2.2. The strong non-Hermiticity limit 

In this limit r € [0, 1) and thus fi^ are kept fixed while the large- limit 
is taken. Therefore both macroscopic and microscopic correlators extend to 
the complex plane. The microscopic rescaling has to be modified, reading 



(23) 



for the eigenvalues and masses. Because //^ does not scale with A^, the 
saddle point analysis for Nf = 2 or more flavors is difficult, as can be seen 
from eq. ()19() . We therefore restrict ourselves to the case Nf = 1, taking 
eq. (EJ as a starting point. If we rescale the mass according to eq. ((231), 
V2Nm = u), and expand for small values of /x^ <C 1, we obtain 



2l,strongip; Oj) = exp[-iV(l + ^^)] exp 



-U! 



1 



lo { V2Nuj{l + -fi' 



(24) 



The evaluation of the second partition function takes a little more care, 
starting from eq. ()13() . Expressing the single Laguerre polynomial as an 
integral, we obtain 



-Z//(r;m) ~ r^exp 



' -Nm?' 


1 


2t 


Jo 



dsexp[-Ns]s^Io { Nm 



-N N 

~ e r exp 



-Nm'' 
2t 



Iq Nm\l - 




(25) 
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In the saddle point expansion is the second step we have aheady used the 
mass rescahng \/TVm = lo, which is appropriate here for (qq)^ = 1/2 with r 
close to unity. If we employ again the relation ()1U() . = 1 — 2fj,'^, insert it 
here and expand in /i^ <^ 1 we arrive at 

Zn,strong{p; i^) = exp[-7V(l + I?)] exp 

x/o [V2Nuj{l + (26) 

We have used that — > exp[— A^/i^], omitting terms of order 0(^^). We 
observe that the leading order (suppression) term in /i^ as well as the term Iq 
exponentially growing in a; match to the given order 0(/i^). We recall that 
the matching of the two macroscopic densities also included an expansion 
of eq. (HJ in to obtain eq. The finite term exp[— a;2/2] only matches 
to zero-th order. 

We find this mapping of the two partitions functions at strong non- 
Hermiticity quite remarkable and can show it to hold also for Nj > 2 and 
small //^ <C 1 Another reason for that is that the strong limit can be 
obtained from the weak limit in taking a ^ oo. While the two partition 
functions may disagree for larger values of /x^, as the expansions indicates, 
this does not exclude a universal behavior of the corresponding microscopic 
correlation functions. A mapping of these quantities may though include a 
nontrivial, /i- and r-dependent unfolding procedure, respectively. 



3. Comparison to lattice data 

In this section we briefly recall the matrix model predictions |lUj and 
their confirmation by quenched QCD lattice data jl2j . A more detailed 
analysis of the data will be published elsewhere 

The quenched microscopic spectral density at weak non-Hermiticity in 
the sector of topological charge v is given by 



/ dte-"'*J,(Vte)^.(Vtr), (27) 

JO 

where the eigenvalues are rescaled according to eq. (IT^ . The density is 
normalized to unity for large arguments on the real line. At strong non- 
Hermiticity we obtain 



Pweak (0 



erf(a) 



\i\ exp 



— (9mO' 



Pstrong (0 



27r 



1 -r^ 



ICI exp 



(28) 
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with the different rescaling eq. (|2,''{|) . Higher order correlation functions can 
be found in jTU]. 

The data we compare with were generated on a 6^ lattice at gauge 
coupling (3 = 5.0. Since staggered fermions at strong coupling away from 
the continuum limit are topology blind as explained in jl9j, we have set 
= in all comparisons below. The two values chosen for the chemical 
potential are /i = 0.006 and /i = 0.2. They correspond to weak and strong 
non-Hermiticity, respectively. We display the same figures as in ref. ^2] 
but with more statistics (17.000 and 20.000 configurations, respectively). 
For simulations with other lattice sizes, confirming the scaling prediction 
eq. (|14|) in the weak limit, we refer again to ^3]. The data are rescaled 



. 5 




MO 



Fig. 2. Cut of the Dirac eigenvalue density along the real axis for /i — 0.006. 



p(i) 



r 



1 



9m (,^) 



Fig. 3. Cuts of the same density along the imaginary axis, with fixed real part of 
the eigenvalues at the first maximum (left) and minimum (right) in fig. [3 

with the mean level spacing d ~ YjV determined as follows. For \i = 0.006 
the data almost lie on the real axis and are therefore ordered according to 
their real part. We thus simply count the number of eigenvalues in a (real) 
interval of given length to determine the average spacing. The same factor 
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d is used to rescale /x^ from eq. (fTl|) and we obtain a = 0.19. In fig. El 
andEleq. (|^7j) is plotted with this value of a versus the rescaled data, and 
we find an excellent agreement. We note that there are no free parameters 
and that no fit has been made^. Fig. is very reminiscent to comparisons 
for // = (see with the difference that we need at least 10 times more 
statistics. This is due to the spread of the data into the complex plane seen 
in fig. EJ which is thus truly testing the value of a. Here, eq. (P7)) is plotted 
for fixed real part 5Re(^) on the first maximum and minimum (saddle) next 
to the origin. Comparing to the histograms that contain the maximum 
or minimum we find an excellent agreement. In jl3j the same picture is 
confirmed for different lattice sizes keeping the value of a fixed. 

p{i) 



9m(^)| 




m 



|Ke(0 



Fig. 4. Cuts of the Dirac eigenvalue density for fi — 0.2 along the imaginary (left) 
and real axis (right). 



At strong non-Hermiticity ^ does not scale with the volume and we 
thus have r = \/l — 2^"^ ~ 0.96 for /i = 0.2. The data are rescaled with the 
inverse square root of the mean level spacing determined as follows. For this 
value of chemical potential the spreading of data into the complex plane is 
such that an ordering according to the real part does no longer make sense. 
For a given set of eigenvalues in our window close to the origin we thus 
determine for each eigenvalue its closest neighbor in the complex plane and 
then average over these distances. After averaging over all configurations 
we obtain the mean level spacing. In fig. 0]left (right) the density eq. (|28)) is 
plotted for for fixed Ke(^) = (9m(^) = 0) and compared to the histograms 
along the imaginary (real) axis. We find very good agreement with the 
parameter free prediction^. The statistical fluctuations are still larger than 
in the weak limit, because for about the same number of configurations 
more eigenvalues have moved away from the real axis. Along the imaginary 
axis the density drops to zero already at ~ ±1.1 in given units, indicating 
that only part of the left histogram is truly in the microscopic regime. 

^ The correction of eq. (IlUt explains the small discrepancy found in |I2| . 
We note that compared to |12| we have refined the determination of the level spacing. 
In addition eq. 1281 has been rescaled, taking properly into account (gg)^ = 1/2. 
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4. Conclusions 

To summarize we liave sliown tlie equivalence of two different matrix 
model realizations 9, 1^ for QCD with chemical potential. The equivalence 
holds for weak non-Hermiticity with several flavors and for strong non- 
Hermiticity with at least one flavor for small chemical potential. It shows 
that the results of can be directly related to QCD. To illustrate this 
we have compared to quenched QCD lattice data for the microscopic Dirac 
operator spectrum and found agreement in both limits, for weak and strong 
non-Hermiticity. 

We close with some remarks on the quenched approximation. While 
the origin of its failure to describe QCD with non-vanishing chemical po- 
tential was clarified in the matrix model is valid beyond the quenched 
approximation as well. It should be mentioned that also for zero chemical 
potential the quenched approximation is not without inconsistencies j2Uj . 
Nevertheless, matrix models do correctly describe unquenched simulations, 
as can be seen e.g. for two colors with dynamical fermions For complex 
spectra so far only quenched or phase quenched predictions with massless 
flavors are available ^Uj. In proving the equivalence we have learned how 
to introduce mass terms. This indicates that with a chemical potential the 
duality between massless flavors and topology proved in 22 may no longer 
hold in general. This seems quite plausible since for example the quenched 
theory with nonzero topology 7^ should not carry a phase, while the 
unquenched theory with the same number of massless flavors Nf and zero 
topology should carry a phase. 
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